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Abstract. We describe a code designed for hot media (T 
> a few 10 4 K), optically thick to Compton scattering. It 
computes the structure of a plane-parallel slab of gas in 
thermal and ionization equilibrium, illuminated on one or 
on both sides by a given spectrum. Contrary to the other 
photoionization codes, it solves the transfer of the con- 
tinuum and of the lines in a two stream approximation, 
without using the local escape probability formalism to 
approximate the line transfer. We stress the importance 
of taking into account the returning flux even for small 
column densities (10 22 cm -2 ), and we show that the es- 
cape probability approximation can lead to strong errors 
in the thermal and ionization structure, as well as in the 
emitted spectrum, for a Thomson thickness larger than 
a few tenths. The transfer code is coupled with a Monte 
Carlo code which allows to take into account Compton 
and inverse Compton diffusions, and to compute the spec- 
trum emitted up to MeV energies, in any geometry. Com- 
parisons with CLOUDY show that it gives similar results 
for small column densities. Several applications are men- 
tioned. 

Key words: Radiative transfer, Galaxies:active, Ultravi- 
olet: galaxies, X-rays:galaxies 



1. Introduction. 

A number of codes were built in the past to compute the 
structure and the emission of photoionized media. With 
time, these codes became more and more sophisticated, 
and able to treat a larger number of situations. In the six- 
ties they were specially designed for planetary nebulae and 
HII regions, i.e. for dilute and optically thin media (except 
in the Lyman continuum) ionized by the thermal radiation 
of a hot star, and were rapidly used also for ionization by a 
non thermal continua extending in the X-ray range, i.e. for 
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the Narrow Line Region of Active Galactic Nuclei (AGN). 
At the same time another generation of codes was devel- 
oped for optically thin hot collisionally ionized media, and 
used for the solar corona, for supernova remnants, and for 
the hot intergalactic gas. At the end of the seventies, ob- 
servations of compact X-ray sources implying the presence 
of a Thomson thick medium incited Ross (1979) to develop 
a radiative transfer method for the continuum using a 
modified version of the Kompaneets equation. In the same 
paper he introduced the so-called "escape probability" ap- 
proximation to take into account diffusion and absorption 
into the lines. This approximation was thereafter amply 
used for the Broad Line Region: this was the beginning of 
the local escape probability era, with the fiducial paper of 
Kwan & Krolik (1981) aimed at studying the Broad Line 
Region of AGN (BLR), immediately followed by a simi- 
lar computation for the atmospheres of X-ray binary stars 
(Kallman and McCray 1982) from which XSTAR derived. 
The most popular of these codes is CLOUDY, designed by 
Fcrland (cf. for instance Ferland & Rees 1988) and contin- 
uously updated since this time (cf. Ferland 1996, Ferland 
et al. 1998). These codes are not only used for relatively 
"cold" clouds like the BLR, but also for warm photoion- 
ized media (Broad Absorption Line region in quasars and 
"Warm Absorbers" in AGN, shock heated media, etc.). 

At the other extreme (i.e. for dense, non irradiated, 
and semi-infinite media) model atmosphere codes have 
been developed since the fifties with an emphasis put onto 
the computation of line transfer, which was treated com- 
pletely, but with the LTE approximation. Non LTE effects 
were introduced in the sixties, and extensively studied in 
the seventies, with the state of the art being beautifully 
exposed in Mihalas book (1978). Since this time many im- 
provments have been made in the numerical methods, in 
particular aiming at taking into account a large number 
of atoms and levels, cf. for introductory reviews Rutten 
1995, Hubeny 1997. 

In the nineties appeared the urgent need for "interme- 
diate" codes, valid for thick or semi-infinite dense media, 
eventually hot, irradiated by a non thermal continuum 
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extending in the hard X-rays, to cover the whole range 
of situations encompassed in AGN and in binary stars. 
Several codes were then built, mainly to compute the 
em-ission spectrum of accretion disks, irradiated or not 
by an X-ray continuum. For irradiated disks, they were 
cither of "photoionization type", using the Kompaneets 
equation (Ross & Fabian, 1993, based on the Ross 1979 
code), or coupling an existing photoionization code to a 
Monte Carlo computation (Zycki et al. 1994), to take into 
account Compton diffusions. On the other hand, sophis- 
ticated model atmosphere codes were transformed to be 
applied to accretion disk structure and emission (Hubeny, 
1990, Hubeny & Hubeny 1997 and 1998), but without the 
external irradiation by an X-ray continuum. 

Owing to the high optical thickness of the medium 
in several frequency ranges, such codes require that the 
transfer of both the continuum and the lines be solved 
in an "exact" way i.e. avoiding approximations such as 
local escape approximation (for the lines) or one stream 
approximation (for the continuum). Since the medium is 
generally dense and sometimes close to LTE, they require 
that all processes and inverse processes be carefully han- 
dled. Being irradiated by an X-ray continuum, the medium 
contains a large number of ionic species, from low to high 
ionization, which should all be introduced in the compu- 
tation. Finally, the medium being hot and thick, not only 
Thomson, but also Compton scattering, should be taken 
into account. 

We have undertaken to build a code in order to sat- 
isfy these requirements. Precisely we have built several 
interconnected codes, which allow more flexibility. The en- 
semble is far from being perfect and still contains several 
approximations which restrict its use, but we intend to 
improve it in the future. 

We have presently four codes. One of them (titan) is 
designed to study the structure of a warm or hot thick 
photoionized gas, and to compute its emission - reflection 
- transmission spectrum from the infrared up to about 
20 KeV. It solves the energy balance, the ionization and 
the statistical equilibria, the transfer equations, in a plane- 
parallel geometry for the lines and continuum. Then, given 
the thermal and ionization stratification, the computation 
of the emitted spectrum from 1 KeV to a few hundre-ds 
KeV is performed with NOAR which uses a Monte-Carlo 
method taking into account direct and inverse Compton 
scattering (it allows also to study various geometries), pe- 
GAS adresses the case of more diluted and thin media, and 
is similar to cloudy with however some differences, and 
IRIS is specially devoted to the computation of the line 
fluxes, including very weak lines, using the most recent 
available atomic data and a full treatment of the statistical 
equilibrium equations for a great number of levels. These 
codes can be used to model a wide variety of astrophysical 
media, particularly those in which the energy transport is 
purely radiative, but also plasmas heated by other mecha- 
nisms (viscosity, shocks, energetic particules...). They have 



already been used in several published papers (Collin- 
Souffrin et al. 1996, Czerny & Dumont 1998, Porquet et 
al. 1998, Abrassart 1999). 

The three codes NOAR, PEGAS and IRIS, are described 
in other papers (Abrassart 2000, Dumont & Porquet 2000). 
In this paper we present titan and we discuss a few simple 
cases, schematized by a slab of gas irradiated on one side 
or on both sides. As far as possible we will perf-orm com- 
parisons with cloudy in the range of parameters where 
both codes can be used, to assess the validity of titan. 
However detailed comparisons are not possible, as they 
use quite different transfer methods in particular. We will 
also show some examples of the capabilities of TITAN in 
the parameter range which is not accessible to CLOUDY. 

titan is mainly designed to determine the structure 
(temperature and ionization state) and the continuum em- 
ission spectrum of a thick hot photoionized slab of gas. 
Owing to the representation of the ions made in the code 
it does not give an accurate determination of the weak 
lines. For a detailed determination of the line spectrum of 
some ions, it should be complemented by IRIS. 

We briefly summarize below the physical processes (Sect. 

2) , the transfer method and the iteration procedure (Sect. 

3) , the thermal equilibrium, (Sect. 4), focussing only on 
the aspects which are not treated in a standard way. The 
coupling of titan and NOAR is briefly described in Sect. 
5, and some applications are presented in Sect. 6. 

2. Physical processes 

In TITAN the physical state of the gas (temperature, ion 
abundances and level populations of all ionic species) is 
computed at each depth, assuming stationary state, i.e. 
local balance between ionizations and recombinations of 
ions, excitations and deexcitations, local energy balance 
(equality of heating due to absorption and cooling due to 
local emission), and finally total energy balance (equality 
between inward and outward fluxes). 

Due to the large range of density and temperature in- 
side the medium, many physical processes play a role at 
some place and should therefore be taken into account. 
Ionization equilibrium equations include radiative ioniza- 
tions by continuum and line photons, collisional ioniza- 
tions and recombinations, radiative and dielectronic re- 
combinations, charge transfer by H and He atoms, the 
Auger effects, and ionizations by high energy electrons 
arising from ionizations by X-ray photons. Energy balance 
equations include free-free, free-bound and line cooling, 
and Compton heating/cooling. 

The emission-absorption mechanisms for the contin- 
uum include free- free and free-bound processes, two-photon 
process, and Thomson scattering. Special care is given to 
recombinations to ground state, which are very important 
from an observational point of view in the X-ray range. 
They are treated differently according to the relative val- 
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ties of kT to the photon energy bin, in order to get an 
accurate frequency dependence. 

Hydrogen and hydrogen-like ions are treated as 6-level 
atoms. Levels 2s and 2p are treated separately, while full 
1-mixing is assumed for higher levels. All processes includ- 
ing collisional and radiative ionizations and recombina- 
tions are taken into account for each level (cf. Mihalas 
1978). Recombinations onto levels n > 5 are not taken 
into account, which amounts to assuming that the higher 
levels are in LTE with the continuum, which is generally 
true in the conditions for which this code is presently used. 
In the future we plane to add several other levels, and to 
sum the contributions of the higher levels as it is done for 
instance in cloudy. Level populations are then obtained 
as usual by matrix inversion. 

In order to save computation time, non H-like ions are 
presently treated with a rough approximation: interlock- 
ing between excited levels is neglected and populations of 
the excited levels are computed separately using a two- 
level approximation. This approximation does not predict 
correctly the details of the line spectrum, since it neglects 
subordinate lines. Nevertheless recombinations onto ex- 
cited states are taken into account in the ionization equi- 
librium and the transfer of these photons is treated in an 
approximative way as proposed by Canfield & Ricchiazzi 
(1980). We assume also that each recombination produces 
a resonant photon after cascades. 

The gas composition include 10 elements (H, He, C, 
N, O, Ne, Mg, Si, S, Fe), and all their ionic species are 
taken into account. Photoionization cross sections are fit- 
ted from Rcilman & Manson (1979) and Band L.M. et 
al. (1990), these values being correct as far as neutral 
and once ionized ions are not concerned, and it is the 
case in hot media. For total radiative and dielectronic 
recombination rates, we use Aldrovandi & Pequignot's 
data (1973). When possible, collisional excitation rates are 
taken from the Daresbury Report (1985). Most of data 
for iron come from Arnaud & Raymond (1992), Kaastra 
& Mewe (1993), and from Fuhr et al. (1988). Ionizations 
by high energy electrons arising from ionizations by X- 
range photons are taken from Bergeron & Souffrin (1973). 
Inverse processes (except dielectronic recombinations) are 
computed through the equations of detailed balance. In 
the case of a gas close to LTE, we neglect dielectronic re- 
combinations for consistency to insure the balance for each 
process. All induced processes are taken into account. 

The equations are not recalled as they have been given 
in previously quoted papers. 

Iron K lines 

These lines require special attention as they are intense 
in Seyfert nuclei, and they will be observed in detail in the 
future with Chandra and XMM. Though the iron K lines 
constitute a complex system described in Band D.L. et 
al. (1990), we assume presently only one "mean" line per 
ion with an oscillator strength equal to 0.4, as suggested 



by Band. More detailed computations are not required as 
far as the Doppler broadening of the lines (or Compton 
broadening as well, see Abrassart 2000) is much larger 
than the distance between the lines, as it is the case in 
AGN. 

In other computations of the line fluxes it is generally 
assumed that resonant trapping of Ka photons of FeXVII 
to Fe XXIII (Ross & Fabian 1993) suppress completely 
these lines when the Thomson thickness of the emitting 
medium is larger than a given value, of the order of 0.02 
(Zycki & Czerny 1994). Here the transfer of these lines is 
handled in a standard way, with an additional term in- 
cluded in the statistical and ionization equations to take 
into account the competition between radiative deexcita- 
tion and the Auger process. The population of the 
upper level of a transition Ka of the ion i is given by: 

N*A kl = ^{Ni-xKi-x + N % B tk J lk ) (1) 

while the ionization equation for the same ion i writes: 

Ni[(Pi + Ki) + B ik J ik {l - (2) 
= Ni+iOi - Ni-iKi-i(l - yi-i). 

A k i and Bi k are the Einstein coefficients of the Ka line, 
Jik is the mean intensity integrated over the line profile, iji 
is the fluorescent yield, a, is the recombination coefficient, 
and Pi (respt. Kj) is the photoionization rate of the ground 
state (respt. of the K-shell) of the ion i. The second term 
on the left side of Eq. || is generally only of the order of 
10% of the total photoionization rate. 

3. Radiation transfer 

An important quantity will be used all along this paper, 
the illumination (or "ionization") parameter. Among sev- 
eral definitions used in the literature, we adopt the follow- 
ing: 



where F mc is the frequency integrated flux incident on one 
side of the slab and «h is the hydrogen number density at 
this surface. Note that it does not preclude the possibility 
of having also a flux incident on the back side of the slab. 
We call the attention on the fact that all along this paper 
we will mimic a semi-isotropic and not a mono-directional 
illumination. 

3.1. Present method for the transfer of the continuum 

We want to study, on one hand media with an optical 
thickness larger than unity in a range of frequency rich in 
emission or in absorption lines, on the other hand media 
with a very inhomogeneous structure. For instance if the 
illuminated side has a temperature of a few 10 6 K, while 
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the temperature is only a few 10 4 K on the back side of 
the cloud, the spectrum emitted by the illuminated side 
is mainly formed in hot layers where the absorption coef- 
ficient is weak at all wavelengths and Thomson diffusion 
is dominant (we shall call it the "reflected spectrum" 
though it does not correspond to pure reflection). 
On the contrary the spectrum emitted by the backside of 
the slab is formed in cold layers where absorption must be 
carefully taken into account (we shall call it the "outward 
spectrum"). Note also that, even if the incident illumi- 
nation is mono-directional, the "transmitted spectrum" 
is almost isotropic since the Thomson thickness, tt, is 
larger than unity, ft is therefore not differentiated from 
the "outward spectrum" even when the incident source is 
not located on the line of sight. It is possible to estimate 
an approximate transmitted spectrum using the expres- 
sion F* ransm = Fl nc exp{-T„), where F™ c is the incident 
flux, in an aim of comparison. 

The transfer is treated with the Eddington two-stream 
approximation, i.e. the intensity is assumed to be constant 
in each hemisphere. The transfer equations can then be 
written : 



J_dll 

y/3 dz 

-l dl~ 



(4) 



V3 dz ■ v v ■ 2' 

where z is the distance to the illuminated edge, k v is the 
absorption coefficient, a is the diffusion coefficient - here 
it is due to Thomson scattering - and e„ is the emissiv- 
ity (all these quantities are local and depend on the fre- 
quency v except a) . Note that this approximation is closer 
to the semi-isotropic case than to the normal case. We 
have already mentioned that it is appropriate for the re- 
flected and for the transmitted flux when the Thomson 
thickness is larger than unity. It is also appropriate for a 
semi-isotropic illumination, if the source of radiation is ex- 
tended (like for instance in the disk-corona model of AGN, 
or in the blob model of Collin-Souffrin et al. 1996). 

The mean intensity J„ is equal to (1+ + 1~)/2, while 
the flux F defined as usual by J I v cos 9 duj is equal to 
(1+ — 7,7) 2it/\/S ; so that the reflected flux is equal to 
J"(Q) 2?r/V3, and the outward flux to I+(H) 2n/V3. The 
optical depth and the total optical depth are defined as : 



T V { Z ) 



/ V3(k,v + <r)dz' and T v = T U {H) 
Jo 



(5) 



where H is the total thickness of the slab, and the source 
function as: 

(e„ + erj„) 



S 



(6) 



Both sides of the slab can be illuminated by an external 
radiation, so the boundary conditions are, at z = 0: 

V3 



and at z = H: 

I-(H) = ^F^ (8) 

where F ^ ack is the back-side flux (it is equal to zero in 
many of the following computations). 

The formal solution of the transfer equations between 
z — 5z and z gives: 



I+(z) = I+(z - Sz)e~ 5r " + e- T " / S u (t)e +t dt. (9) 

J t u — 8t v 

Assuming that S v (t) is nearly constant in the interval 
z, z + Sz , one gets: 

I+(z) = I+(z-6z)e- ST » + 1 ~ S 2 [S v {z-6z)+S v {z)]{VS) 

with a similar equation for I~(z). 

Guessing the initial value of 1^(0) fails because calcu- 
lations diverge, unless this initial value is provided with an 
extreme precision: the predicted value should differ from 
the real one at all frequencies by less than 10~ 6 . This 
problem has also been pointed out by Coleman (1993). 
We therefore start to compute I„(z) from the back side 
with: 

I-{z) = I-(z+Sz)e- Sr -+ 1 ~ 6 2 ^ \S v {z+5z)+S v {z)}{11) 

o Thus we calculate S u (z) and the mean intensity J u (z) 
from Eq. || and: 



T ( x _ Ij-(z-Sz) + I-{z + 6z) „ STv 

•Jy\Z) — „ e 



(12) 



(1-e-^) 



[S v {z - Sz) + 2S v (z) +S v (z + Sz)} 



4 + (0) 



2tt 



TTiinc 



(7) 



assuming a constant Sz (actually it is not constant so the 
formulae are more complicated but do not deserve to be 
given here). 

The procedure is then the following. We divide the sl- 
ab into a set of plane-parallel layers with the density tih 
given in each layer. Note that the code can easily be cou- 
pled with a prescription for the density or for the pressure, 
such as a constant pressure, or a pressure given by hydro- 
static equilibrium. Presently the slab is divided in about 
300 layers. The geometrical thickness of each layer is sm- 
aller close to the surfaces than in the middle of the slab 
because of the rapid variation of the physical paramet- 
ers (the optical thickness, the temperature, the ionization 
state). 

- For each layer, starting from the illuminated side, where 
ljj~(0) is given by Eq. [?], the ionization and thermal balance 
equations are solved. The source function, the opacity and 
the emissivity, are computed. As they depend on J Vl we 
use Eq. |l3| with the value of I^(z — Sz) of the previous 
layer and the value of I~ (z + Sz) provided by the previous 
iteration and given by Eq. |ll] (see below). 
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Fig. 1. The reflected continuum vF v , in erg cm~ 2 s -1 , 
after a number of iterations varying from 100 to 600, dis- 
played all 100 iterations, for the reference model. We see 
that the spectrum does not vary at all after 400 iterations. 

- I„(z) is transferred through the layer according to Eq. 
|To| , while S v (z) and t v {z) are buffered for each layer and 
each frequency. 

- when the back side of the cloud is reached, new values 
of Iy{z) are calculated, starting from the back side where 
I~ (H) is given by Eq. pi and using Eq. |ll| with the previ- 
ous values of S v {z) and r v {z) ; actually we use the values 
given by several previous iterations to accelerate the pro- 
cedure. 

- the whole calculation is repeated until convergence. It 
is stopped when energy balance is achieved for the whole 
slab (i.e. when the flux entering on both sides of the slab 
is equal to the flux coming out from both sides). 

The drawback of this method is that it requires a large 
number of iterations to get a complete convergence of the 
model. The iteration starts by assuming given values for 
l~(z). If this initialisation is not adequate, convergence 
problems are encountered. We have chosen for the first 
iteration I~(z) equal to I„(z) multiplied by a factor of 
the order of unity, which works reasonably well. Actually 
the convergence of ly(z) (as well as of the temperature 
and the ionization state) is reached within less than 1% 
after only a relatively small number of iterations in the 
layers corresponding to about one Thomson thickness. 
This property is interesting as far as the outward spec- 
trum is negligible, which is the case in the EUV and X- 
ray range. It can also be used to determine the ionization 
state and the temperature of the layers where Compton 
diffusion takes place (cf. Sect. 5 about coupling with the 
code noar). 

Fig. ^ displays the reflected continuum obtained after 
different numbers of iterations (for clarity we prefer to 
show here the continuum and not the complete spectrum) . 
The model is a slab with the following characteristics: 

- a constant density 10 12 cm~ 3 , 

- a column density 10 26 cm" 2 ; 
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Fig. 2. Variation of the surface temperature, in K, with 
the column density, for a slab illuminated on one side by 
a power law continuum with £ = 100 erg cm s . The 
surface temperature computed by CLOUDY for £ multiplied 
by vo is shown for comparison. 

- it is illuminated on one side by an incident power law 
continuum proportional to v^ 1 from 0.1 eV to 100 KeV, 

- the illumination parameter is £ — 10 3 erg cm s , 

- there is no illumination on the other side; 

- the abundances are (in number): H: 1, He: 0.085, C: 3.3 
10" 4 , N: 9.1 10" 5 , O: 6.6 10" 4 , Ne: 8.3 10" 5 , Mg: 2.6 10" 5 , 
Si: 3.3 10~ 5 , S: 1.6 10~ 5 , Fe: 3.2 lO" 5 . 

In the following we will call this our "reference model" . 
All along the paper we will use the same density and the 
same spectral distribution of the illuminating radiation, so 
we will not specify these parameters anymore. Note that, 
for this peculiar spectral distribution, £ and the ratio of 
incident ionizing photon number to hydrogen densities, U, 
are related by : £ = 8.21 U ln(i/ max / V m ; n ) where u max and 
^min are the maximum and minimum photon frequencies. 

We want to stress here the importance of the back- 
ward intensity. It cannot be neglected as it participates to 
the ionization, to the level population, and to the energy 
balance. For instance in a purely scattering medium the 
returning flux varies as Fi nc j^—^ ; so when tt is much 
larger than unity, all the incident flux is reflected. Thus 
the temperature of the first layer at the illuminated side 
(which depends on the sum of incident plus returning radi- 
ation) should increase with the column density for a given 
incident flux (except possibly in cases of very high values 
of the illumination parameter). 

Fig. H shows the dependence of the surface tempera- 
ture on the column density N. The model is similar to 
the reference one except that the illumination parameter 
is equal to 100 erg cm s _1 and the column density is vary- 
ing. We see that the surface temperature is constant up to 
about N — 10 22 cm~ 2 , but after it increases from 2 10 5 K 
to 2.7 10 5 K for N increasing from 10 22 to 10 24 cm" 2 . It 
saturates when a large fraction of the slab becomes close 
to LTE, and the returning flux becomes constant. With 
CLOUDY (version 9004) the surface temperature does not 
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Fig. 3. Dependence of the surface temperature, in K, on 
the ratio /~(0)// + (0) imposed a priori. 

depend on the column density (for a correct comparison 
the illumination parameter in CLOUDY must be multiplied 
by a factor V3 to account for the semi-isotropic approxi- 
mation used in titan). 

Fig. H displays the dependence of the surface temper- 
ature on an assigned ratio /~(0)/J + (0) independent of 
the frequency. We see that a modest returning flux (20%) 
is sufficient to account for the 60% increase of the surface 
temperature reached for TV = 10 24 cm -2 , with respect to a 
very small column density. To understand this behaviour, 
one should remember that the spectral distribution of the 
reflected flux contains a much larger proportion of EUV 
radiation than the incident continuum, in particular in 
the lines. So the illumination parameter is increased by a 
larger amount than the ratio of the integrated reflected 
over incident continuum. We can also guess that the in- 
crease of the temperature at the illuminated surface with 
respect to cloudy will be compensated by a decrease at 
the back side, because of energy conservation (cf. Fig. |l4j ). 

3. 2. Present method for line transfer 

Radiative transfer of the lines is treated in the same way 
as the continuum, the line profile being represented by 
several frequencies distributed symmetrically around the 
line center. 

We assume that even if large macroscopic velocities are 
present, they do not play any role in the line transfer, in 
other words that the thickness of the slab is much smaller 
than the scale length of the velocity gradient. We also 
assume that the lines do not overlap, which is valid as far 
as multiplets are treated globally, and if the microscopic 
turbulent velocity is not much larger than the thermal 
velocity (doing this we neglect line fluorescence due to the 
overlap of two lines, like the Hell Balmer/3 line and the 
hydrogen Lya line, whose wavelengths differ only by a 
fraction 4 10~ 4 ). 

Finally we assume presently complete redistribution in 
all the lines. We intend to include partial redistribution for 
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Fig. 4. Profile of I v and /+ for OVIII Lya, for the refer- 
ence model. The curves are labelled by the value of z/H. 



resonance lines like Lya in the near future. For the mo- 
ment partial redistribution in these lines can be mimicked 
by assuming complete redistribution in a pure Doppler 
profile, as the line photons diffuse more in space than in 
frequency, owing to the absence of enlargment of the lower 
level. A Doppler profile can also be used for strongly in- 
terlocked lines, and in any case, it does not differ strongly 
from a Voigt profile when To is smaller than 1/a so the 
medium is optically thin in the damping wings. Note that 
partial redistribution is less important in dense media, like 
irradiated accretion disks in AGN, than in dilute media, 
owing to collisional redistribution. 

The line profile </>„ is a symmetrical Voigt function. 
The Doppler width Az^d and the damping constant depend 
on the temperature and on the density, and consequently 
vary with z. It may happen that the temperature decreases 
by two orders of magnitude from the illuminated to the 
back side, so Afd decreases by one order of magnitude. 
This must be taken into account to choose the points in 
the profile. To solve the statistical equilibrium equations 
one needs to know the total line intensity weighted by the 
profile J J v cf)(v)dv, and to treat the thermal and ionization 
equilibria one simply integrates the line intensity over the 
frequencies J J v dv. For intense lines the first integral is 
dominated by the Doppler core, and the second by the 
wings. In the case of a saturated line with strong wings, 
the set of frequencies chosen to describe the line must 
cover a range from a few tenths A^j to a few hundreds 
A^d and at least 10 points are necessary for a correct 
representation of the profile. 
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Fig. 5. Variation of p, P e (ro)/2 and /3 e ft as functions of 
the distance to the illuminated side for several intense 
lines, for the reference model. 



Fig. |i| displays for the reference model the variation 
with z/H of the OVIII Lya profile in both directions, I~ 
and J+, J v being the half sum of these profiles. We see 
that the profile varies considerably from the surface to 
the deeper layers. 

Now we are able to show the difference between the 
transfer treatment and the escape probability approxima- 
tion for the lines discussed in the Appendix. In this aim, we 
have compared the divergence flux p ( Eq. [j"7| ) with the es- 
cape probability towards the illuminated side P e (ro), and 
with an approximate expression proposed by Rees et al. 
(1989) to take into account continuum absorption in the 
statistical and in the energy balance equations : 



'cff 



+ 



Kl 



[Pe(T ) + Pe(T - T ) 



Kl 



Kl 



(13) 



where f3 c s is the total effective escape probability, k\ and 
Kjq are the absorption coefficients respectively at the line 
center and in the underlying continuum, To is the total 
optical depth of the slab at the line center, and P e {To) 
(respt. P e (To — To)) are the escape probabilities towards 
the illuminated (respt. the back) side. 



Lya 
Ha 

Hell Lya 
Hell 1640 A 
CIV 1549 A 
OIII 306 A 
OIV 280 A 
OVI 1034 A 
OVII 22 A 
OVm Lya 
OVm Ly(3 
Fe K XXIII 
Fe K XXIV 
FeXXV 6.7 keV 



-f .275 A 
-4.44 A 
9.47 A 
1.76 A 
-4.15 A 
0.042 A 

1.07 A 
1.83 A 
0.12 A 

2.08 A 
-0.0012 A 

-32 eV 
303 eV 
170 eV 



5.61 A 
21.4 A 
6.3 10 4 A 
10.4 A 
-11.6 A 
1230 A 
569 A 
-1.61 A 



Table 1. Equivalent widths of a few intense lines (nega- 
tive for absorption, positive for emission) in the reflected 
and in the outward spectra, for the reference model. 



Pe(T~o)/2 and (3 c g have been computed with titan and 
compared to p in Fig. ^, as functions of the distance z to 
the illuminated surface, for the reference model. 

We see on this figure that (3 e g and P e (ro)/2 are almost 
always identical. It is first due to the fact that the slab 
is very thick, so P c (to) is much larger than P e (To — tq) 
except very close to the back side, and second because k c 
is always <C ki, as all the lines displayed on the figure are 
intense. But /3 c h differs considerably from p in the region 
where the lines are formed (i.e. the region where p is not 
equal to zero and where the ion is present) . In particular p 
is negative in an important fraction of the line formation 
region for Ha, La, and OVIII La, while /3 e ff is a positive 
quantity. 

As a consequence the line spectrum is not correctly 
computed when escape probabilities are used. This can be 
seen on Table [l] which gives the equivalent widths (nega- 
tive if the line is in absorption, positive if it is in emission) 
of some intense lines from both sides of the slab. Several of 
these lines are in absorption, either from the illuminated 
or from the dark side, which is impossible with the escape 
probability approximation. 

Presently the radiative transfer is solved for 300 fre- 
quencies in the continuum from 0.01 eV to 25 keV, in- 
cluding all the ionization edges, and for 440 spectral lines 
including the fluorescence Fe K lines, the integrals over the 
line profiles being achieved using 15-point Gauss-Legendre 
quadrature, together with a change of variable. 

Remarks on the convergence of the method 

It is well-known that the lambda iteration method used 
here (actually mainly for historical reasons, because we st- 
arted to build a nebular code with the line escape probabil- 
ity approximation) converges extremely slowly for intense 
lines in stellar atmospheres, and we intend to replace it 
by the Accelerated Lambda Iteration method as soon as 
possible. The convergence of J u is achieved after 600 iter- 
ations, for all lines except a few ones and for all continua. 
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Fig. 6. Profiles of some intense lines emitted by the bright 
side after a number of iterations varying from 100 to 600, 
displayed all 100 iterations, for the reference model. 

This can be checked on Fig. ^ which displays for the same 
lines as Fig. ^| the profiles of the lines emitted by the bright 
side, after different number of iterations, for the reference 
model. As a consequence, the structure of the slab (T, 
ionization) does not change at all after a few hundreds 
iterations. Indeed our case is more favorable than stellar 
atmospheres. The explanation is double: 1. the relatively 
small density of the medium compared to a stellar atmo- 
sphere, which decreases comparatively the influence of col- 
lisional excitations, and 2. the fact that in a photoionized 
plasma the temperature is low for a high ionization degree, 
so the excited levels lie at very high energies compared to 
the thermal energy. Consequently the populations of these 
highly ionized ions are smaller than at Boltzman equilib- 
rium with the ground level, and the lack of convergence 
for excited level populations does not have a strong influ- 
ence on the ionization state of the corresponding ion. We 
have also checked that the same structure is obtained for 
different initial boundary conditions. All this implies that 
the continuum emitted from both sides (including for in- 
stance the Lyman discontinuity) is completely converged, 
and only the few non converged lines are not well com- 
puted, in particular those emitted by the back side (but 
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Fig. 7. Dependence of the local cooling, heating, and net 
cooling-heating rates, in erg cm 3 s _1 , and of the temper- 
ature, in K, on the distance to the illuminated side z, for 
the reference model. 

they are negligible in the whole energy balance). Since 
we do not claim to compute a correct line spectrum be- 
fore the subordinate lines are implemented in the code, we 
think that this approximation is sufficient for our present 
purpose. 

4. Energy balance 

titan has been designed for stationary radiatively heated 
media, so the energy conservation implies the balance of 
radiative flux. However it can be used also for other cases 
out of radiative equilibrium, for instance by imposing the 
temperature of the medium, or by adding a non radiative 
energy in each layer, such as the viscous flux of an accre- 
tion disk. We describe here how the code works in a purely 
radiative case. 

We first stress the fact that the energy balance should 
be realized both locally and globally. 

The local energy balance equation writes: 
r dp 

J -^-dv = =n e n H [-T tot + A tot ], (14) 

where n e n#r tot and n e riH A tot are the total heating and 
cooling rates per unit volume. The temperature of a layer 
is computed through an iteration process until r to t and 
Atot are equal. 

As an example, Fig. [?] gives the dependence of the 
cooling and heating rates, and of the equilibrium tem- 
perature, on the distance from the illuminated side z, for 
the reference model. We note that in this model Compton 
heating-cooling is not very important, and that bound- free 
contributes as a heating term, compensated by a cooling 
term due to lines and free-free processes. Note also that 
a quasi LTE is reached at a few Thomson depths in this 
model (it is reached at a smaller optical depth for a sm- 
aller value of £). Above this value of z the temperature 
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Fig. 8. Dependence of the total line and continuum fluxes and of the fluxes of a few intense lines on the number of 
iterations, for the reference model. 



is almost constant, except close to the back side where it 
decreases more rapidly due to the leakage of photons. 

The global energy balance is reached through an it- 
eration procedure, when the total flux emitted by both 
sides of the slab, F t o ™, is equal to the total flux incident 
on both sides, F^fr. Actually this is the most important 
convergence criterium and it requires many iterations. To 
illustrate the problem, Fig. || displays the flux of a few 
intense emission lines and of the total flux in the lines and 
continuum, for the reference model, as functions of the 
number of iterations. We see that after about 200 itera- 
tions the global energy balance is reached within about 
one percent. 

In the same way, Fig.|| shows that the energy balance 
with titan is reached within one percent after 200 it- 
erations for two models with a column density equal to 
10 25 cm -2 , and with an illumination parameter equal to 
£ = 100 and £ = 1000 erg cm s -1 . CLOUDY or xstar have 
been used several times in this range of parameters (Mar- 
tocchia & Matt, 1996, Zycki et al. 1994, for instance), but 
it is not sure that the global energy balance is achieved 
owing to the escape treatment of the lines, as discussed 
above. 

Finally we must recall that the neglect of subordinate 
lines (except in hydrogen-like ions) could lead to a sub- 
stantial increase of the energy losses and consequently to 
a change of the thermal balance, as these lines have a 
smaller optical thickness than resonant lines, and there- 
fore can escape more easily from the deep layers. Since in 
a photoionized plasma the temperature is relatively low 
and the excited level populations of highly ionized ions 
are generally small, one could expect only a small influ- 
ence on the overall spectrum. 

To check this point we compare in Figs. [l(] and [ll] 
the temperature and the emitted spectrum for the refer- 
ence model, with and without taking into account inter- 
locking between excited levels in hydrogen like ions (note 
that interlocking is included for hydrogen and helium in 
both cases). These ions are the most abundant ones in 
the region emitting the reflected spectrum, in our refer- 
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Fig. 10. Temperature in K versus z with (full line) and 
without (dashed lines) full five-level interlocked hydrogen- 
like ions, for the reference model. 

ence model. Fig. [n] shows that the temperature is slightly 
shifted in layers corresponding to a Thomson thickness 
of a few. As a consequence the emitted spectrum is also 
different, but not by an important amount. We can infer 
a similar behaviour for smaller illumination parameters, 
even if the abundant ions are less ionized, as the tempera- 
ture is also lower. For a higher illumination parameter the 
lines contribute less to the energy balance. 

Since the detailed computed line spectrum of hydrogen- 
like ions is different with the two treatments, contrary to 
the overall spectrum, it is mandatory to treat correctly the 
excited levels for all the ions in order to get a correct line 
spectrum. Note that in this respect CLOUDY is presently 
doing much better than titan, as it includes subordinate 
lines and more levels for each ion. 

Compton heating and cooling 

The Compton heating-cooling term is equal to: 

r - A = ^-—4k(T Comp - T) [ 4nJ u dv (15) 

where Tc om p is the Compton temperature. We give here 
this well-known equation only to recall that it depends (as 
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Fig. 11. Reflected and outward spectra, with (full line), and without (dashed lines) full five- level interlocked hydrogen- 
like ions, for the reference model. The spectra vF v are displayed with a line resolution of 30, in erg cm~ 2 s _1 . 



well as Tcomp) on the spectral distribution of the local flux, 
so it requires a careful computation of this flux (including 
the lines). 

Below 25 KeV, J u is provided directly by titan in 
each layer. Above 25 KeV, J„ is not calculated by titan, 
but by the code NOAR, as the energy shift of the photons 
due to Compton scattering is no more negligible. In this 
case the results from NOAR are also used to compute the 
exact value of the Compton heating-cooling term. This is 
discussed in more detail in the next section. 



5. Coupling with the Monte Carlo code 

In order to compute the effect of Compton scattering on 
the high energy part of the reprocessed radiation, and par- 
ticularly on the reflected spectrum above 1 keV, and to 
determine precisely the Compton heating-cooling rate en- 
tering in the energy balance equation, we need another nu- 
merical approach. The transfer scheme in titan is based 
on frequency arrays, and inelastic scattering would in- 
troduce a redistribution on energy bins making the code 
much heavier. In this aim Abrassart (2000) developped a 



Monte Carlo code, NOAR. The asset of such an approach 
is also that it enables to investigate an arbitrary geometry 
and to determine the angular dependance of the observed 
spectra. These aspects were difficult to tackle with a code 
such as titan. Moreover, it allows to easily extract time 
variability information. 

NOAR is described in detail in Abrassart (2000). We ex- 
plain here only how it is coupled with titan. It shares the 
same opacity data and includes the same set of 102 ions, 
for the 10 elements of highest cosmic abundances. Given 
all the ionisation degrees provided by titan, it computes 
absorption cross sections in each layer. Free-free absorp- 
tion is also taken into account. NOAR includes the recom- 
bination continua of hydrogen and helium like ions, but 
does not include line emission, except fluorescence lines. 
The proper yields to account for the competition with the 
Auger effect are used for these lines. Fluorescence of Fe 
XVII- XXIII is suppressed by resonnant trapping. For the 
pseudo fluorescence of H and He-like species, a simplify- 
ing prescription is adopted which includes only the two 
most probable outcomes of a K shell photoionization, i.e. 
direct recombination on ground level or L-K transition. 
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The spectral distribution of the recombination continuum 
is determined by the local temperature in the "current" 
slab. 

noar takes into account direct and inverse Compton 
scattering. The method used for modelling Compton scat- 
tering is basically the one described in Pozdniakov, Sobol 
& Sunyaev (1983) and in Gorecki & Wilczewski (1984), 
with a different use of statistical weight, because of the 
competition with photoelectric absorption. When the high 
energy cut-off of the incident continuum is of the order of 
100 KeV, the reflected spectrum above 10 KeV exhibits a 
"Compton hump" peaking at about 30 KeV (Lightmann 
& White 1988). This spectral region of optimal continuum 
albedo occurs where the sum of photoelectric absorption 
and Compton losses is minimum. This feature depends 
weakly on the ionization state, although its prominence 
in an observed X-ray spectra, and notably its level with 
respect to the iron line, is sensitive to £ (see Fig. |l9|) as 
the soft X-ray albedo. Note that the shape and extension 
of the hump depends on the high energy cut-off of the 
incident spectrum, which is 100 KeV all along this paper. 

Fluorescent lines are significantly Compton broadened 
for rather high illumination parameter, above £ ~ 10 3 
erg cm s _1 . The broadening is asymmetyric, the profile is 
skewed toward the red. At still higher £ and temperature, 
a blue wing appears, due to the inverse Compton effe- 
ct. The extend and importance of both wings depends 
on the optical depth of the scattering medium. The red 
wing depends only on the frequency of the line (i.e. on the 
ionization state of iron) whereas the blue one only depends 
on the temperature of the scattering medium. 

titan and NOAR lead to almost identical reflected spec- 
tra in the 1-20 KeV range (they differ only by Compton 
broadening of lines and of photoelectric edges), and noar 
allows to obtain the spectrum in the higher energy range. 
This spectrum is used for the global energy balance. 

Another use of NOAR is to provide titan with the local 
Compton gains and losses in each layer. This is necessary, 
because Compton heating-cooling is dominated by energy 
losses of photons > 25 Kev, in particular for high values of 
£ when Compton heating-cooling is dominating the energy 
balance and should be computed in an exact way. 

As an example Fig. [l2| gives the dependence of the ra- 
tio (rcomp - Acomp)/£ on the Thomson depth, for different 
values of £, the other parameters being the same as in the 
reference model. We see that the ratio does not depend 
on the illumination parameter if £ < 10 4 erg cm s — 1 . This 
is expected since the Compton heating-cooling is domi- 
nated by high energy photons, which do not care about 
the ionization state. It also means that inverse Compton 
is negligible as a cooling process (though it has an influ- 
ence on the line spectrum), otherwise ^-j^ would depend 
on the gas temperature and therefore on £. This is due 
to the relatively low temperature of the gas compared to 
the mean energy of the photons. This property breaks for 
very high illumination parameters, since low energy pho- 




Fig. 12. Dependence of the ratio rc ° mp g Ac " mp on the 
Thomson depth, for different values of £ in erg cm s _1 , 
the other parameters being the same as in the reference 
model. 
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Fig. 13. Dependence of the ratio rc ° mp Ac "°p on the 
Thomson depth, for a neutral gas and for an isotropic or 
perpendicular illumination. Except £ the parameters are 
the same as in the reference model. 

tons arc then able to penetrate in the deepest layers and 
to play a role in the cooling rate, and the temperature of 
the gas is close to Tc om p- 

The deposition of energy shown on Fig. [l^ corresp- 
onds to a semi-isotropic distribution of the incident radi- 
ation. It is different in the case of a mono-directional ra- 
diation. This is shown on Fig. [l3| which displays (Tcomp — 
Acomp)/£ for a neutral gas in the case of an isotropic and 
of a perpendicular illumination. When the radiation field 
is semi-istropic there is a rapid decrease of for small 
values of tt , due to the small distance covered by the pho- 
tons. This is compensated by a smaller (by about a factor 
two) deposition of energy in the deeper layers. 

The Compton heating-cooling rate obtained with NOAR 
is then fitted analytically as a function of z, which is trans- 
ferred to TITAN. 
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6. Some applications 



A full grid of models will be published elsewhere (Dumont 
& Abrassart 2000) . Here we discuss only a few interesting 
cases. 



6.1. Optically thin hot medium: the Warm Absorber in 
AGN 

The Warm Absorber (WA) is a hot medium located on 
the line of sight of the X-ray source in AGN. It is respon- 
sible of the absorption edges of OVII and OVIII observed 
in their soft X-ray spectrum in about 50% of them. It 
might also constitute the "mirror" invoked in the Unified 
Scheme to account for the polarized Broad Lines observed 
in Seyfert 2 galaxies (Antonnucci & Miller 1985). It has 
been modelled in many papers, using ion (Netzer 1993, 
1996), xstar (Krolik & Kriss 1995), cloudy (Nicastro, 
Fiore & Matt 1999 for the most recent paper), pegas and 
IRIS (Porquet et al. 1999). This medium has a column 
density of 10 21 to 10 24 cm -2 , and a typical illumination 



parameter of £ = 100 erg cm s _1 . Porquet et al. have 
shown that its density should be at least 10 10 cm -3 to 
avoid the emission of too strong coronal lines. 

As the WA is relatively thin (r„ <C 1 for the continuum 
in a large range of frequencies), it is possible to use the line 
escape probability formalism. This is however limited to 
the case where r„ is smaller than unity in the continuum 
underlying intense lines. This specific problem will be dis- 
cussed in detail in the paper describing the code PEGAS 
(Dumont & Porquet 2000). A correct computation of the 
backward flux is also required. We have shown indeed that 
the temperature depends on the intensity of this backward 
flux. 

We have run several models corresponding to the con- 
ditions of the WA, i.e. for £ = 100 erg cm s _1 and column 
densities varying from 10 22 to 10 24 cm~ 2 (the other pa- 
rameters being the same as in the reference model). We 
have compared the results to those obtained with CLOUDY 
for the same models, to check that they are similar at low 
column densities, and to determine at which column den- 
sities they begin to differ strongly. Note that a detailed 
comparison between titan and CLOUDY is not possible, 
as these codes differ both in the transfer method and in the 
atomic data (for instance owing to the semi-isotropic illu- 
mination, a given column density corresponds to a larger 
optical thickness with titan than with CLOUDY) , so only 
trends can be obtained. 

Fig. [l4| displays the dependence of the temperature on 
z, obtained with titan and with CLOUDY. As discussed 
previously, it is higher with titan at the illuminated side 
for a high column density, owing to the returning flux, and 
it is lower at the dark side. 

The ionization state (of fundamental importance for 
the study of the WA) and the emission spectrum, are also 
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Fig. 15. Fractional abundances of Carbon, Oxygen, and 
Iron, computed by titan (full lines) and CLOUDY (dashed 
lines), as functions of the distance to the illuminated side, 
for £ = 100 erg cm s^ 1 and a column density equal to 10 24 
cm~ 2 . Note that the absissa scale for Iron differs from the 
others. 
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Fig. 16. Variation of the flux at 54.9 cV (the OIII-OIV 
ionization edge) as a function of the depth, and compari- 
son with the transmitted flux. 

different with titan and with cloudy at high column 
densities. Fig. [f5| shows the ionization state of Carbon, 
Oxygen and iron for a column density of 10 24 cm -2 , com- 
puted with TITAN and with CLOUDY. While the ionization 
state close to the illuminated side is quite similar (in par- 
ticular if one remembers that there is a factor V3 in the 
optical depth), it is very different in the deepest layers, 
where CLOUDY leads to a smaller degree of ionization than 
titan (CIII and OIII instead of CIV and OIV). This effect 
cannot be explained by the difference in optical thickness, 
which acts in the opposite direction. It is certainly due to 
the different treatments of the diffuse continuum, and per- 
haps of the lines. With titan we have checked that the 
diffuse continuum at large depths is sufficient to ionize 
OIII into OIV. To illustrate this point Fig. |l6| displays the 
variation of the flux at 54.9 eV (the OIII-OIV ionization 
edge) as a function of the depth. It shows that the dif- 
fuse flux is still very intense and sufficient to ionize OIII 
at large depths (contrary to the transmitted flux which 
is very weak). Note that this result depends strongly on 
the way the spectral distribution of the diffuse flux is com- 
puted. If it would be reduced to a unique frequency at the 
recombination frequency of Hc ++ , 54.4 eV, (as it is gen- 
erally the case with the "On The Spot" approximation), 
there will be no photons able to ionize OIII. 

Fig. [I?] displays the reflected and outward spectra, for 
column densities 10 23 and I0 24 cm~ 2 . We see that the 
overall shape of the continuum is very similar with titan 
and with cloudy. However, some differences in the de- 
tailed features appear for the largest column density. For 
instance the Lyman edge is in emission in the reflected 
spectrum with both codes, but it is larger with titan 
(presumably because the temperature is higher in the em- 
ission layers). In the outward spectrum it is present as a 
very weak absorption with CLOUDY, while it is in emission 
with titan. 

The line spectra are obviously also different. There are 
much more lines with CLOUDY, but on the other hand the 
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Fig. 18. Temperature in K versus z for different values of 
£ in erg cm s _1 , the other parameters being the same as 
in the reference model. 

fewer lines of titan are more intense. This is expected 
since subordinate lines are not taken into account exce- 
pt for hydrogen-like ions in titan. As mentioned before, 
some lines appear in absorption with titan, both in the 
reflected and in the outward spectrum. 

6.2. Hot medium optically thick to Compton scattering: 
the UV-soft X-ray emitting medium in AGN 

TITAN is specially designed for Compton thick photoion- 
ized media like those commonly assumed to emit the UV/soft 
X-ray continuum in AGN and to produce through Comp- 
ton reflection the FeK ~7 KeV line and the 30 KeV hump 
observed in many Seyfert 1 galaxies. 

Although the exact nature of this medium is not known 
(an irradiated accretion disc, as proposed by Ross & Fabian 
1993, Zycki et al. 1994, to quote only the first papers 
on the subject, or a clumpy Compton thick medium, cf. 
Collin-Souffrin et al. 1996), its characteristics and phys- 
ical state are comparable: a shell of gas with a column 
density of at least 10 25 cm" 2 , a temperature of 10 5 -10 6 K 
due to radiative heating, and a density spanning a range 
from 10 12 to 10 15 cm -3 . The corresponding value of £ is 
typically 300 to 3000 erg cm s _1 . Our reference model is 
therefore representative of this medium. In the case of an 
irradiated accretion disk, it has been mimicked as a sh- 
ell of constant density irradiated from above by a power 
law continuum, and from below by a black body radia- 
tion (actually it is not the best way to take into account 
the presence of the underlying viscously heated disk, cf. 
Rozanska et al. 2000). 

Let us consider first a slab illuminated on one side by 
a semi-isotropic radiation. Fig. ^ displays T versus z, for 
different values of £. The other parameters are the same 
as in the reference model. It is interesting to notice that 
the quasi constant temperature regime (corresponding to 
LTE) is reached for increasing values of z when £ increases. 
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Fig. 17. Reflected and outward spectrum, for £ = 100 erg cm s _1 and column densities equal to 10 23 and 10 24 cm~ 2 . 
The results obtained with CLOUDY are shown for comparison in dashed lines. The incident flux of CLOUDY is multiplied 
by a factor y3 to account for the semi-isotropy of the intensity in titan, and the emitted fluxes are multiplied by an 
additional factor 10 for a better visibility. The spectra are displayed with a line resolution of 300. 



For instance it is reached for z ~ 10 13 cm, corresponding 
to 7T ~ 6 for £ = 1000 erg cm s _1 . It means that any com- 
putation aiming at giving the reflected spectrum of such 
a model should solve correctly the thermal and ionization 
equilibrium until these deep layers. 



Fig. 19 shows the reflected spectrum for the same mod- 



els. Note that the resolution adopted here (100) does not 
correspond to a region where the lines are broadened by 
large velocity field, either turbulent or organized (rota- 
tion) . We use this high resolution only to show more clearly 
the details of the spectrum. 

We see that the overall shape of the spectrum com- 
puted with titan and NOAR is very similar in the range 1- 
20 KeV. However the detailed spectral features are differ- 
ent, particularly for the higher values of £. For instance the 
complex feature near 7 KeV (which is made of a mixture 
of several Iron edges and of Iron lines dominated by low 
ionization stages for £ = 100 erg cm s _1 , and by FeXXV 
and FeXXVI for £ = 10 4 ) erg cm s _1 is strongly modified 
by Compton scattering. For £ = 10 4 erg cm s _1 the line 
has a large red wing due to direct Compton scattering, 
and a weak blue wing due to inverse Compton scattering. 
The absorption edge is also completely erased. The effect 
is smaller for £ = 10 3 erg cm s _1 , and almost absent for 
£ = 100 erg cm s _1 , as discussed above. Similar results 
have been obtained in many papers dealing with the X- 



ray spectrum of Seyfert galaxies (cf. for instance Ross et 
al. 1999). 

Let us now compare the results obtained by us with 
TITAN, and by Zycki et al. (1994), for a slab of constant 
density 10 14 cm -3 , illuminated by a power law v~ 9 be- 
tween 10 eV and 100 KeV on one side, and by a black 
body Tbb = 10 5 K on the other side. We will compare 
the results for £(Zycki et al.)=300, which corresponds to 
£(titan)^ 200 erg cm s _1 , owing to our semi-isotropic 
illumination, and to a slightly different definition of £. 

Like us, Zycki et al. used a Monte-Carlo method to 
take into account Compton scattering. The thermal and 
ionization state of the gas is calculated apart with the code 
XSTAR, using like CLOUDY the escape probability approx- 
imation. Fig. ^o] gives T versus tt obtained by us and by 
Zicky et al. The same result is observed as with CLOUDY: 
our temperature is larger near the surface of the slab. It 
reaches Tbb at tt=0.7. For a larger value of £, the connec- 
tion with the underlying black body would be reached in 
even deeper layers, as discussed above. For tt > 0.7 the 
temperature becomes smaller than Tbb, and joins Tbb 
only at tt ~ 20. This behaviour is however dependent on 
the thickness of the slab, it is why this representation of 
the disk is not meaningful (Rozanska et al. 2000). On the 
contrary Zycki et al. impose the gas temperature to be 

equal to min(T B B, ^thermal balance)- 
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Fig. 19. Reflected spectrum for different values of £ in 
erg cm s , the other parameters being the same as in 
the reference model; full lines: results of noar; dashed 
lines: results of titan. The spectra are displayed with a 
resolution of 100. 
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Fig. 20. Temperature in K versus tt for the "accretion 
disk" model (see text), computed by us (full line), and by 
Zycki et al. (dashed line). 
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Fig. 21. Reflected spectrum for the "accretion disk" 
model computed by us (full line), and by Zycki et al. 
(dashed line). Zycki et al. spectrum is shifted for clarity. 
The resolution is 20 for the titan spectrum. 

Fig. |l] displays the reflected spectrum, computed by 
us and by Zycki et al. The resolution is degraded to 20, to 
fit the energy bins of Zycki et al. The overall shape of the 
spectra are roughly similar. The major discrepancy is the 
OVIII edge, in emission with our computation and in ab- 
sorption with Zycki et al.'s one. We think that it is due to 
the higher surface temperature and to the lower temper- 
ature in the deeper layers obtained in our computation. 

7. Conclusion 

We have shown that the coupling of titan and NOAR al- 
lows to compute the structure and the emission of hot 
Compton thick irradiated media in an unprecendented 
way. First it solves consistently both the global and the lo- 
cal energy balance, which is impossible in a thick medium 
with codes handling the line transfer with the escape prob- 
ability approximation, as all present photoionization codes 
do. We have also shown the importance of the returning 
flux (which is neglected in photoionization codes) even 
for relatively low column densities. Second it takes into 
account in an exact way inverse and direct Compton scat- 
tering, both in the energy balance and in the computation 
of the emitted spectrum. Finally, it allows to treat any ge- 
ometry, open or closed. 

Although the problem of the convergence process is not 
as drastic as in stellar atmospheres, since we are able to 
get complete converged structure and spectrum in a still 
reasonable computing time, a most urgent improvement of 
the code is to accelerate the convergence process through 
the use of the Accelerated Lambda Iteration method. This 
will not only allow to get the results in a much smaller 
time, but also to get convergence for the few lines which 
are still not converged after about 10 3 iterations. 

Then the following improvements of TITAN will be to 
take into account subordinate lines in solving a multi-level 
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atom for all ions and to bring up to date the atomic data. 
More detailed line spectra for some abundant ions will 
be obtained through coupling with the code IRIS. The L 
shells of Iron which are already taken into account for 
FcXVII to FcXXII will be implemented as well as a better 
representation of the Iron K lines. A few elements will be 
added to the already ten existing ones. 
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Appendix : Previous methods 

7.1. Transfer of the continuum 

It is not always clear which approximations are used in the 
different codes. Basically the computation of the diffuse 
radiation field seems to be similar in all codes, except in 
Ross & Fabian (1993) and in subsequent works using this 
code, where it is computed using the Kompaneets equa- 
tion, cloudy (Ferland 1996), xstar (Kallman & Krolik, 
1995 for the last version), ion (Netzer 1990, 1993), use 
a modified version of the "on the spot" approximation, 
which amounts to assuming a kind of escape probability 
for the diffuse continuum, and one stream approximation 
in the transfer equation. These approximations are correct 
only for a continuum optical thickness smaller than unity 
(effective or total opacity) , and when the properties of the 
cloud do not vary considerably between the point where 
a photon is emitted and the point where it is reabsorbed. 
In particular one very important requirement for studying 
hot and thick media is to take into account the radiation 
returning from the backside of the slab, even when it is not 
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illuminated. As the medium is optically thick in a large 
frequency range, this radiation is intense and modifies the 
physical state of the whole slab, including the illuminated 
side(see Fig. @). 



7.2. Line transfer: the escape probability 
approximation 

All previous photoionization codes, except that of Collin- 
Souffrin and Dumont (S.) (1986), treat line transfer by the 
so-called local escape probability formalism, which uses 
the probability that a line photon emitted at a given point 
can escape in a single flight from the cloud, assuming 
that the rest of the cloud is homogeneous and has the same 
properties as the emitting layer. It amounts to identifying 
the divergence flux of the statistical equilibrium equations 
with the escape probability intervening in the line emerg- 
ing flux. This approximation is valid only in the case of 
complete frequency redistribution, absence of line inter- 
locking, and if the medium is homogeneous. Its use 
can have severe consequences on the emission line spec- 
trum and on the energy balance when these conditions 
are not fulfilled as it was often discussed, see for example 
by Collin & S.Dumont (1986) or Elitzur (1984). We in- 
troduce here the divergence flux in the aim to show that 
our results are completely different from those obtained 
by the escape probability formalism. 

To simplify the discussion, let us consider in the rest of 
this section a simple two-level atom. The level population 
balance writes: 



n u (A u i + B u i J J v ip v dv + n c C u i) 
= ni(Bi u / J v $ v dv + n c Ci u ) 



(16) 



where <p v and tj) v are the absorption and emission line pro- 
files, A u i, B u i and Bi u , C u i and Ci u are the usual radiative 
(Einstein) and collisional excitation and deexcitation co- 
efficients, and 3 V is the angle averaged intensity. 

Let us now define the divergence flux of the transition, 

pul- 



Pul 



{n u (A u i + B u i J Jv^vdv) - ni(B lu J J v 4> v dv)} 



(17) 



n u A u i 

With this definition Eq. [l7| becomes simply: 

n u (n c C u i + PuiAui) = n t n c Q u . (18) 

Let us assume the line is the only contributor to ra- 
diation. The energy emitted locally in the line in a small 
volume of surface 1 cm~ 2 and length dz writes: 

dF = hv [n u (A u i + B u i J J v i\) v dv) - niB lu J J v <\> v dv\ dz 
= Pui"n u Aul hv dz — p u iAiredz (19) 



where e is the emissivity coefficient n u A u ihu / A-k . In- 
tegrating this equation over depth, one gets the emerging 
line flux: 



F = 4w / pui e dz 



(20) 



According to this equation the divergence flux seems 
to be equal to the probability that, once emitted at a 
distance z from the surface, a photon can escape from the 
medium, which is called the "escape probability" P e . 

From this fact, if one makes the approximation that the 
emission and the absorption profiles are equal, one com- 
putes this escape probability by integrating the attenua- 
tion over the line profile: 

p ^ 1 Jo J if n u A ul (j> v exp{-T§ v l[i)dv dfi 

2 J if n « A ul 4>v dv 

where /i = cos#, and it gives the usual expressions for 
complete redistribution in a Voigt profile: 



P e (ro) ~ max 



1 + 2r 0v /7r/n(ro + 1)' 3 V t ^k 



(22) 



where To is the optical thickness at the line center and a 
is the usual damping constant. This expression can thus 
replace p u i in Eqs. [l8] and |2(]. 

With this method the frequency dependent transfer eq- 
uations are successfully replaced by frequency integrated 
ones, and solved consistently with the statistical equat- 
ions. 

Several estimates of the escape probability are used in 
the literature (see for example Kwan & Krolik, 1981, Rees, 
Netzer & Ferland, 1989, Kallman & McCray, 1982, Ko & 
Kallman, 1994), to account for partial redistribution, for 
line or continuum interlocking, and for continuum absorp- 
tion. For instance several ways to take into account the 
continuum opacity (or the overlapping of two lines) have 
been proposed (Elitzur & Netzer, 1984, Netzer, Elitzur & 
Ferland, 1985). In the case of the Broad Line Region of 
AGN, which is optically thin for the continuum underly- 
ing the main lines, Collin-Souffrin et al. (1981), Avrett 
& Loeser (1988), have discussed the influence of differ- 
ent approximations for partial or complete redistribution, 
and how they compare to an exact transfer treatment, and 
Collin-Souffrin & Dumont S. (1986) have shown that the 
escape probability approximation is roughly valid (it can 
however lead to a decrease by a factor as large as two of 
the ratio of a lower to a higher transition, such as Ha/H/3 
or La/Ha). 

However these approximations are not valid if r cont is 
of the order of unity at the line frequency and if the line 
photons are created in one place and absorbed in another 
place where the physical conditions are different. This hap- 
pens for instance for EUV and soft X-ray lines created in 
the hot region and absorbed in ionizing He + ions in a 
cooler region. 
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If the local escape probability approximation breaks 
down, it has a severe consequence: the emitted line spec- 
trum and the thermal and the ionization balance are not 
correctly computed. For instance the divergence flux takes 
frequently negative values in our computations, even for 
intense lines, and even in the reflected spectrum. This is 
not allowed with escape probabilities. We illustrate this 
discussion with some examples in Sect. 3.2. 



